Confounder boar (Igf & animal welfare)

Author

Anja Eggert & Anne-Marie Galow

Published

March 2, 2025

R Libraries

library(tidyverse)     # tidy universe
library(patchwork)     # combine plots
library(lmerTest)      # mixed model
library(car)           # ANOVA
library(emmeans)       # post hoc
library(performance)   # model performance
set.seed(1989)

Data

Read data

Data processing done in file “igf-biomarker-summary-stats.qmd”.

load("./data/data-processed.RData")
table(dat.w$husbandry, dat.w$boar)
              
               Astein Bestein Beton Capello Faro Flausist Sacher Slavon Slovan
  conventional      1       0     6       6    0        2      2      1      2
  ecological        4       2     9       1    1        1      0      0      0
              
               Tonkar
  conventional      0
  ecological        2

Statistical model: general design

With a simple model structure we solely test for a confounder effect of the boars. We include a random animal effect and fit a linear mixed model with the lmerTest package in R.

contr = lmerControl(optimizer   = "bobyqa",
                    optCtrl     = list(maxfun = 10000000),
                    calc.derivs = FALSE)

Birth weight

Model

hist(dat.w$mean.bodyweight,
     breaks = 30)

mod.bodyweight <- lmerTest::lmer(mean.bodyweight ~ 
                                   boar +
                                   # random intercept for sows
                                   (1 | sow),
                                 data    = dat.w,
                                 REML    = TRUE,
                                 control = contr)
summary(mod.bodyweight)
Linear mixed model fit by REML. t-tests use Satterthwaite's method [
lmerModLmerTest]
Formula: mean.bodyweight ~ boar + (1 | sow)
   Data: dat.w
Control: contr

REML criterion at convergence: 10.2

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.23173 -0.57132  0.08495  0.58636  1.93133 

Random effects:
 Groups   Name        Variance Std.Dev.
 sow      (Intercept) 0.003039 0.05513 
 Residual             0.055767 0.23615 
Number of obs: 40, groups:  sow, 15

Fixed effects:
             Estimate Std. Error       df t value Pr(>|t|)    
(Intercept)   1.28340    0.10830 29.99958  11.850 7.63e-13 ***
boarBestein   0.19258    0.20261 29.99985   0.951    0.349    
boarBeton     0.12817    0.12416 27.88475   1.032    0.311    
boarCapello  -0.05304    0.14080 27.95471  -0.377    0.709    
boarFaro     -0.16233    0.26503 29.87555  -0.612    0.545    
boarFlausist  0.04580    0.17568 28.20866   0.261    0.796    
boarSacher    0.28907    0.20241 29.86622   1.428    0.164    
boarSlavon    0.08854    0.26514 29.96194   0.334    0.741    
boarSlovan   -0.11074    0.20261 29.99998  -0.547    0.589    
boarTonkar    0.04312    0.19963 23.10452   0.216    0.831    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) brBstn borBtn brCpll boarFr brFlss brSchr boarSlavn
boarBestein -0.534                                                    
boarBeton   -0.860  0.467                                             
boarCapello -0.752  0.402  0.654                                      
boarFaro    -0.408  0.218  0.363  0.307                               
boarFlausst -0.606  0.324  0.527  0.462  0.248                        
boarSacher  -0.534  0.286  0.467  0.417  0.218  0.324                 
boarSlavon  -0.408  0.218  0.357  0.307  0.167  0.248  0.218          
boarSlovan  -0.534  0.285  0.464  0.402  0.218  0.338  0.285  0.245   
boarTonkar  -0.515  0.276  0.451  0.388  0.211  0.312  0.276  0.211   
            boarSlovn
boarBestein          
boarBeton            
boarCapello          
boarFaro             
boarFlausst          
boarSacher           
boarSlavon           
boarSlovan           
boarTonkar   0.275   
round(drop1(mod.bodyweight, test = 'Chisq'), 3)
Single term deletions using Satterthwaite's method:

Model:
mean.bodyweight ~ boar + (1 | sow)
     Sum Sq Mean Sq NumDF  DenDF F value Pr(>F)
boar  0.409   0.045     9 28.805   0.815  0.606
car::Anova(mod.bodyweight,
           test.statistic = "Chisq",
           type = 2)
Analysis of Deviance Table (Type II Wald chisquare tests)

Response: mean.bodyweight
      Chisq Df Pr(>Chisq)
boar 7.3391  9     0.6019
car::Anova(mod.bodyweight,
           test.statistic = "F",
           type = 2)
Analysis of Deviance Table (Type II Wald F tests with Kenward-Roger df)

Response: mean.bodyweight
          F Df Df.res Pr(>F)
boar 0.7149  9 26.965 0.6907

Model diagnostics

performance::check_model(mod.bodyweight)

Plot

plot.bodyweight <- dat.w  |> 
  # make plot
  mutate(jit = jitter(as.numeric(boar), 0.3)) |>  
  ggplot(aes(y   = mean.bodyweight)) +
  geom_boxplot(aes(x   = boar, 
                   fill = husbandry),
               col = "black",
               outlier.shape = NA, width = 0.5) +
  geom_jitter(aes(x   = boar,
                  fill = husbandry),
              col = "black",
              position = position_jitterdodge(jitter.width  = 0.15,
                                              dodge.width   = 0.5),
              size = 2, alpha = 0.7,
              show.legend = FALSE) +
  scale_fill_manual(labels = c("Conventional",
                               "Ecological"),
                    values = c("#FFA040", "#008000")) +
  scale_y_continuous(lim = c(0, 2), 
                     breaks = seq(0, 2, 0.5),
                     minor_breaks = seq(0, 2, by = 0.1) ) +
  labs(x = "",
       y = "Mean birth weight [kg]",
       fill = "Husbandry") +
  my_theme +
  theme(legend.position = "top") +
  guides(y = guide_axis(minor.ticks = TRUE),
         fill = guide_legend(override.aes = list(linetype = 0, size  = 5)))
plot.bodyweight

Total piglets

Model

hist(dat.w$total.piglets,
     breaks = 30)

mod.piglets <- lmerTest::lmer(total.piglets ~ 
                                   boar +
                                   # random intercept for sows
                                   (1 | sow),
                                 data    = dat.w,
                                 REML    = TRUE,
                                 control = contr)
summary(mod.piglets)
Linear mixed model fit by REML. t-tests use Satterthwaite's method [
lmerModLmerTest]
Formula: total.piglets ~ boar + (1 | sow)
   Data: dat.w
Control: contr

REML criterion at convergence: 173

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.2833 -0.4447  0.0667  0.5248  1.3179 

Random effects:
 Groups   Name        Variance Std.Dev.
 sow      (Intercept) 5.536    2.353   
 Residual             9.320    3.053   
Number of obs: 40, groups:  sow, 15

Fixed effects:
             Estimate Std. Error      df t value Pr(>|t|)    
(Intercept)   17.3468     1.6209 29.6138  10.702 1.08e-11 ***
boarBestein   -3.4548     3.0186 28.3752  -1.145    0.262    
boarBeton     -0.2354     1.7414 16.6375  -0.135    0.894    
boarCapello   -1.2092     1.9754 16.5531  -0.612    0.549    
boarFaro       5.3426     3.7983 18.4212   1.407    0.176    
boarFlausist  -1.8109     2.4825 18.2299  -0.729    0.475    
boarSacher    -1.5936     2.8986 18.2757  -0.550    0.589    
boarSlavon     3.4623     3.8503 21.4644   0.899    0.379    
boarSlovan     3.2054     2.9944 25.6667   1.070    0.294    
boarTonkar     0.5989     2.6953 12.5549   0.222    0.828    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) brBstn borBtn brCpll boarFr brFlss brSchr boarSlavn
boarBestein -0.516                                                    
boarBeton   -0.809  0.477                                             
boarCapello -0.676  0.361  0.616                                      
boarFaro    -0.396  0.219  0.407  0.285                               
boarFlausst -0.562  0.300  0.511  0.431  0.237                        
boarSacher  -0.513  0.276  0.484  0.458  0.220  0.307                 
boarSlavon  -0.400  0.216  0.376  0.282  0.172  0.256  0.216          
boarSlovan  -0.514  0.273  0.462  0.363  0.215  0.401  0.274  0.367   
boarTonkar  -0.417  0.226  0.398  0.294  0.180  0.244  0.226  0.176   
            boarSlovn
boarBestein          
boarBeton            
boarCapello          
boarFaro             
boarFlausst          
boarSacher           
boarSlavon           
boarSlovan           
boarTonkar   0.222   
round(drop1(mod.piglets, test = 'Chisq'), 3)
Single term deletions using Satterthwaite's method:

Model:
total.piglets ~ boar + (1 | sow)
     Sum Sq Mean Sq NumDF  DenDF F value Pr(>F)
boar 78.053   8.673     9 18.337    0.93  0.522
car::Anova(mod.piglets,
           test.statistic = "Chisq",
           type = 2)
Analysis of Deviance Table (Type II Wald chisquare tests)

Response: total.piglets
      Chisq Df Pr(>Chisq)
boar 8.3744  9     0.4969
car::Anova(mod.piglets,
           test.statistic = "F",
           type = 2)
Analysis of Deviance Table (Type II Wald F tests with Kenward-Roger df)

Response: total.piglets
        F Df Df.res Pr(>F)
boar 0.84  9 23.759 0.5877

Model diagnostics

performance::check_model(mod.piglets)

Plot

plot.piglets <- dat.w  |> 
  # make plot
  mutate(jit = jitter(as.numeric(boar), 0.3)) |>  
  ggplot(aes(y   = total.piglets)) +
  geom_boxplot(aes(x   = boar, 
                   fill = husbandry),
               col = "black",
               outlier.shape = NA, width = 0.5) +
  geom_jitter(aes(x   = boar,
                  fill = husbandry),
              col = "black",
              position = position_jitterdodge(jitter.width  = 0.15,
                                              dodge.width   = 0.5),
              size = 2, alpha = 0.7,
              show.legend = FALSE) +
  scale_fill_manual(labels = c("Conventional",
                               "Ecological"),
                    values = c("#FFA040", "#008000")) +
  scale_y_continuous(lim = c(0, 22), 
                     breaks = seq(0, 22, 5),
                     minor_breaks = seq(0, 22, by = 1) ) +
  labs(x = "",
       y = "Number of piglets",
       fill = "Husbandry") +
  my_theme +
  theme(legend.position = "top") +
  guides(y = guide_axis(minor.ticks = TRUE),
         fill = guide_legend(override.aes = list(linetype = 0, size  = 5)))
plot.piglets

Proportion males

Model

hist(dat.w$prop.males,
     breaks = 30)

mod.males <- lmerTest::lmer(prop.males ~ 
                                   boar +
                                   # random intercept for sows
                                   (1 | sow),
                                 data    = dat.w,
                                 REML    = TRUE,
                                 control = contr)
boundary (singular) fit: see help('isSingular')
summary(mod.males)
Linear mixed model fit by REML. t-tests use Satterthwaite's method [
lmerModLmerTest]
Formula: prop.males ~ boar + (1 | sow)
   Data: dat.w
Control: contr

REML criterion at convergence: 246.1

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.3908 -0.4465  0.0000  0.5513  2.2803 

Random effects:
 Groups   Name        Variance Std.Dev.
 sow      (Intercept)   0.0     0.00   
 Residual             152.7    12.36   
Number of obs: 40, groups:  sow, 15

Fixed effects:
             Estimate Std. Error       df t value Pr(>|t|)    
(Intercept)   55.3165     5.5268  30.0000  10.009 4.48e-11 ***
boarBestein    0.4527    10.3397  30.0000   0.044    0.965    
boarBeton    -10.7701     6.3818  30.0000  -1.688    0.102    
boarCapello   -8.9729     7.2363  30.0000  -1.240    0.225    
boarFaro     -17.2213    13.5379  30.0000  -1.272    0.213    
boarFlausist  -1.2389     9.0253  30.0000  -0.137    0.892    
boarSacher   -10.0959    10.3397  30.0000  -0.976    0.337    
boarSlavon     2.5782    13.5379  30.0000   0.190    0.850    
boarSlovan     0.9993    10.3397  30.0000   0.097    0.924    
boarTonkar     1.5548    10.3397  30.0000   0.150    0.881    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) brBstn borBtn brCpll boarFr brFlss brSchr boarSlavn
boarBestein -0.535                                                    
boarBeton   -0.866  0.463                                             
boarCapello -0.764  0.408  0.661                                      
boarFaro    -0.408  0.218  0.354  0.312                               
boarFlausst -0.612  0.327  0.530  0.468  0.250                        
boarSacher  -0.535  0.286  0.463  0.408  0.218  0.327                 
boarSlavon  -0.408  0.218  0.354  0.312  0.167  0.250  0.218          
boarSlovan  -0.535  0.286  0.463  0.408  0.218  0.327  0.286  0.218   
boarTonkar  -0.535  0.286  0.463  0.408  0.218  0.327  0.286  0.218   
            boarSlovn
boarBestein          
boarBeton            
boarCapello          
boarFaro             
boarFlausst          
boarSacher           
boarSlavon           
boarSlovan           
boarTonkar   0.286   
optimizer (bobyqa) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')
round(drop1(mod.males, test = 'Chisq'), 3)
Single term deletions using Satterthwaite's method:

Model:
prop.males ~ boar + (1 | sow)
     Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
boar 1172.5  130.28     9    30   0.853  0.575
car::Anova(mod.males,
           test.statistic = "Chisq",
           type = 2)
Analysis of Deviance Table (Type II Wald chisquare tests)

Response: prop.males
     Chisq Df Pr(>Chisq)
boar 7.677  9      0.567
car::Anova(mod.males,
           test.statistic = "F",
           type = 2)
Analysis of Deviance Table (Type II Wald F tests with Kenward-Roger df)

Response: prop.males
         F Df Df.res Pr(>F)
boar 0.733  9 26.984 0.6756

Model diagnostics

performance::check_model(mod.males)

Plot

plot.males <- dat.w  |> 
  # make plot
  mutate(jit = jitter(as.numeric(boar), 0.3)) |>  
  ggplot(aes(y   = prop.males)) +
  geom_boxplot(aes(x   = boar, 
                   fill = husbandry),
               col = "black",
               outlier.shape = NA, width = 0.5) +
  geom_jitter(aes(x   = boar,
                  fill = husbandry),
              col = "black",
              position = position_jitterdodge(jitter.width  = 0.15,
                                              dodge.width   = 0.5),
              size = 2, alpha = 0.7,
              show.legend = FALSE) +
  scale_fill_manual(labels = c("Conventional",
                               "Ecological"),
                    values = c("#FFA040", "#008000")) +
  scale_y_continuous(lim = c(0, 80), 
                     breaks = seq(0, 80, 20),
                     minor_breaks = seq(0, 80, by = 5) ) +
  labs(x = "",
       y = "Proportion males [%]",
       fill = "Husbandry") +
  my_theme +
  theme(legend.position = "top") +
  guides(y = guide_axis(minor.ticks = TRUE),
         fill = guide_legend(override.aes = list(linetype = 0, size  = 5)))
plot.males

Proportion stillborn

Model

hist(dat.w$prop.stillborn,
     breaks = 30)

hist(log(dat.w$prop.stillborn+1),
     breaks = 30)

mod.stillborn <- lmerTest::lmer(log(dat.w$prop.stillborn+1) ~ 
                                   boar +
                                   # random intercept for sows
                                   (1 | sow),
                                 data    = dat.w,
                                 REML    = TRUE,
                                 control = contr)
summary(mod.stillborn)
Linear mixed model fit by REML. t-tests use Satterthwaite's method [
lmerModLmerTest]
Formula: log(dat.w$prop.stillborn + 1) ~ boar + (1 | sow)
   Data: dat.w
Control: contr

REML criterion at convergence: 106.9

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-1.32881 -0.42267  0.00086  0.47963  1.82883 

Random effects:
 Groups   Name        Variance Std.Dev.
 sow      (Intercept) 0.426    0.6527  
 Residual             1.138    1.0668  
Number of obs: 40, groups:  sow, 15

Fixed effects:
             Estimate Std. Error       df t value Pr(>|t|)  
(Intercept)   1.25893    0.54075 29.86540   2.328   0.0269 *
boarBestein  -1.11985    1.01007 29.69135  -1.109   0.2765  
boarBeton     0.31743    0.59564 24.42302   0.533   0.5989  
boarCapello  -0.06806    0.67573 24.42948  -0.101   0.9206  
boarFaro      1.72993    1.29344 26.18582   1.337   0.1926  
boarFlausist  1.25723    0.84656 25.36699   1.485   0.1498  
boarSacher   -0.20312    0.98734 26.09434  -0.206   0.8386  
boarSlavon   -0.79010    1.30408 27.67283  -0.606   0.5495  
boarSlovan    1.14679    1.00678 29.18801   1.139   0.2639  
boarTonkar   -1.30818    0.93170 20.48662  -1.404   0.1753  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) brBstn borBtn brCpll boarFr brFlss brSchr boarSlavn
boarBestein -0.525                                                    
boarBeton   -0.828  0.477                                             
boarCapello -0.702  0.374  0.626                                      
boarFaro    -0.402  0.219  0.395  0.291                               
boarFlausst -0.577  0.308  0.516  0.440  0.240                        
boarSacher  -0.523  0.281  0.480  0.447  0.220  0.312                 
boarSlavon  -0.405  0.217  0.371  0.290  0.170  0.250  0.217          
boarSlovan  -0.524  0.279  0.463  0.375  0.217  0.382  0.279  0.335   
boarTonkar  -0.445  0.239  0.412  0.319  0.188  0.262  0.239  0.185   
            boarSlovn
boarBestein          
boarBeton            
boarCapello          
boarFaro             
boarFlausst          
boarSacher           
boarSlavon           
boarSlovan           
boarTonkar   0.237   
round(drop1(mod.stillborn, test = 'Chisq'), 3)
Single term deletions using Satterthwaite's method:

Model:
log(dat.w$prop.stillborn + 1) ~ boar + (1 | sow)
     Sum Sq Mean Sq NumDF  DenDF F value Pr(>F)
boar 14.102   1.567     9 25.813   1.377  0.249
car::Anova(mod.stillborn,
           test.statistic = "Chisq",
           type = 2)
Analysis of Deviance Table (Type II Wald chisquare tests)

Response: log(dat.w$prop.stillborn + 1)
      Chisq Df Pr(>Chisq)
boar 12.391  9     0.1922
car::Anova(mod.stillborn,
           test.statistic = "F",
           type = 2)
Analysis of Deviance Table (Type II Wald F tests with Kenward-Roger df)

Response: log(dat.w$prop.stillborn + 1)
          F Df Df.res Pr(>F)
boar 1.2425  9 25.027 0.3149

Model diagnostics

performance::check_model(mod.stillborn)

Plot

plot.stillborn <- dat.w  |> 
  # make plot
  mutate(jit = jitter(as.numeric(boar), 0.3)) |>  
  ggplot(aes(y   = prop.stillborn)) +
  geom_boxplot(aes(x   = boar, 
                   fill = husbandry),
               col = "black",
               outlier.shape = NA, width = 0.5) +
  geom_jitter(aes(x   = boar,
                  fill = husbandry),
              col = "black",
              position = position_jitterdodge(jitter.width  = 0.15,
                                              dodge.width   = 0.5),
              size = 2, alpha = 0.7,
              show.legend = FALSE) +
  scale_fill_manual(labels = c("Conventional",
                               "Ecological"),
                    values = c("#FFA040", "#008000")) +
  scale_y_continuous(lim = c(0, 40), 
                     breaks = seq(0, 40, 10),
                     minor_breaks = seq(0, 40, by = 2) ) +
  labs(x = "",
       y = "Proportion stillborn [%]",
       fill = "Husbandry") +
  my_theme +
  theme(legend.position = "top") +
  guides(y = guide_axis(minor.ticks = TRUE),
         fill = guide_legend(override.aes = list(linetype = 0, size  = 5)))
plot.stillborn

Proportion spreizer

Model

hist(dat.w$prop.spreizer,
     breaks = 30)

hist(log(dat.w$prop.spreizer+1),
     breaks = 30)

mod.spreizer <- lmerTest::lmer(log(dat.w$prop.spreizer+1) ~ 
                                   boar +
                                   # random intercept for sows
                                   (1 | sow),
                                 data    = dat.w,
                                 REML    = TRUE,
                                 control = contr)
boundary (singular) fit: see help('isSingular')
summary(mod.spreizer)
Linear mixed model fit by REML. t-tests use Satterthwaite's method [
lmerModLmerTest]
Formula: log(dat.w$prop.spreizer + 1) ~ boar + (1 | sow)
   Data: dat.w
Control: contr

REML criterion at convergence: 79.6

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.3649 -0.1482 -0.1482  0.0000  2.4683 

Random effects:
 Groups   Name        Variance Std.Dev.
 sow      (Intercept) 0.0000   0.0000  
 Residual             0.5938   0.7706  
Number of obs: 40, groups:  sow, 15

Fixed effects:
               Estimate Std. Error         df t value Pr(>|t|)    
(Intercept)   9.692e-16  3.446e-01  3.000e+01   0.000 1.000000    
boarBestein  -1.128e-15  6.447e-01  3.000e+01   0.000 1.000000    
boarBeton     1.142e-01  3.979e-01  3.000e+01   0.287 0.776102    
boarCapello   1.052e+00  4.512e-01  3.000e+01   2.331 0.026660 *  
boarFaro      1.751e+00  8.441e-01  3.000e+01   2.075 0.046699 *  
boarFlausist  6.603e-01  5.628e-01  3.000e+01   1.173 0.249878    
boarSacher   -9.460e-16  6.447e-01  3.000e+01   0.000 1.000000    
boarSlavon   -9.362e-16  8.441e-01  3.000e+01   0.000 1.000000    
boarSlovan    2.609e+00  6.447e-01  3.000e+01   4.046 0.000336 ***
boarTonkar    9.173e-01  6.447e-01  3.000e+01   1.423 0.165100    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) brBstn borBtn brCpll boarFr brFlss brSchr boarSlavn
boarBestein -0.535                                                    
boarBeton   -0.866  0.463                                             
boarCapello -0.764  0.408  0.661                                      
boarFaro    -0.408  0.218  0.354  0.312                               
boarFlausst -0.612  0.327  0.530  0.468  0.250                        
boarSacher  -0.535  0.286  0.463  0.408  0.218  0.327                 
boarSlavon  -0.408  0.218  0.354  0.312  0.167  0.250  0.218          
boarSlovan  -0.535  0.286  0.463  0.408  0.218  0.327  0.286  0.218   
boarTonkar  -0.535  0.286  0.463  0.408  0.218  0.327  0.286  0.218   
            boarSlovn
boarBestein          
boarBeton            
boarCapello          
boarFaro             
boarFlausst          
boarSacher           
boarSlavon           
boarSlovan           
boarTonkar   0.286   
optimizer (bobyqa) convergence code: 0 (OK)
boundary (singular) fit: see help('isSingular')
round(drop1(mod.spreizer, test = 'Chisq'), 3)
Single term deletions using Satterthwaite's method:

Model:
log(dat.w$prop.spreizer + 1) ~ boar + (1 | sow)
     Sum Sq Mean Sq NumDF DenDF F value Pr(>F)   
boar 17.747   1.972     9    30   3.321  0.006 **
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
car::Anova(mod.spreizer,
           test.statistic = "Chisq",
           type = 2)
Analysis of Deviance Table (Type II Wald chisquare tests)

Response: log(dat.w$prop.spreizer + 1)
      Chisq Df Pr(>Chisq)    
boar 29.886  9  0.0004587 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
car::Anova(mod.spreizer,
           test.statistic = "F",
           type = 2)
Analysis of Deviance Table (Type II Wald F tests with Kenward-Roger df)

Response: log(dat.w$prop.spreizer + 1)
          F Df Df.res  Pr(>F)  
boar 3.0219  9 26.984 0.01249 *
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Model diagnostics

performance::check_model(mod.spreizer)

Emmeans & Effect sizes

Emmeans:

emm <- emmeans(mod.spreizer,
        pairwise ~ boar, 
        data    = dat.w, 
        adjust  = "tukey",
        lmer.df = "satterthwaite",
        type    = "response")

emm
$emmeans
 boar     response    SE df lower.CL upper.CL
 Astein      0.000 0.345 30   -0.505    1.021
 Bestein     0.000 0.545 30   -0.671    2.043
 Beton       0.121 0.223 30   -0.253    0.683
 Capello     1.863 0.834 30    0.579    4.189
 Faro        4.762 4.440 30    0.194   26.800
 Flausist    0.935 0.861 30   -0.220    3.802
 Sacher      0.000 0.545 30   -0.671    2.043
 Slavon      0.000 0.771 30   -0.793    3.825
 Slovan     12.580 7.400 30    3.463   40.324
 Tonkar      1.503 1.360 30   -0.178    6.615

Degrees-of-freedom method: satterthwaite 
Confidence level used: 0.95 
Intervals are back-transformed from the log(mu + 1) scale 

$contrasts
 contrast           estimate    SE df t.ratio p.value
 Astein - Bestein      0.000 0.645 30   0.000  1.0000
 Astein - Beton       -0.114 0.398 30  -0.287  1.0000
 Astein - Capello     -1.052 0.451 30  -2.331  0.3995
 Astein - Faro        -1.751 0.844 30  -2.075  0.5576
 Astein - Flausist    -0.660 0.563 30  -1.173  0.9711
 Astein - Sacher       0.000 0.645 30   0.000  1.0000
 Astein - Slavon       0.000 0.844 30   0.000  1.0000
 Astein - Slovan      -2.609 0.645 30  -4.046  0.0106
 Astein - Tonkar      -0.917 0.645 30  -1.423  0.9103
 Bestein - Beton      -0.114 0.580 30  -0.197  1.0000
 Bestein - Capello    -1.052 0.618 30  -1.702  0.7856
 Bestein - Faro       -1.751 0.944 30  -1.856  0.6966
 Bestein - Flausist   -0.660 0.703 30  -0.939  0.9936
 Bestein - Sacher      0.000 0.771 30   0.000  1.0000
 Bestein - Slavon      0.000 0.944 30   0.000  1.0000
 Bestein - Slovan     -2.609 0.771 30  -3.385  0.0531
 Bestein - Tonkar     -0.917 0.771 30  -1.190  0.9683
 Beton - Capello      -0.938 0.353 30  -2.658  0.2357
 Beton - Faro         -1.637 0.796 30  -2.057  0.5689
 Beton - Flausist     -0.546 0.487 30  -1.121  0.9785
 Beton - Sacher        0.114 0.580 30   0.197  1.0000
 Beton - Slavon        0.114 0.796 30   0.143  1.0000
 Beton - Slovan       -2.494 0.580 30  -4.300  0.0055
 Beton - Tonkar       -0.803 0.580 30  -1.385  0.9227
 Capello - Faro       -0.699 0.824 30  -0.849  0.9969
 Capello - Flausist    0.391 0.532 30   0.736  0.9990
 Capello - Sacher      1.052 0.618 30   1.702  0.7856
 Capello - Slavon      1.052 0.824 30   1.277  0.9514
 Capello - Slovan     -1.557 0.618 30  -2.520  0.2985
 Capello - Tonkar      0.134 0.618 30   0.218  1.0000
 Faro - Flausist       1.091 0.890 30   1.226  0.9619
 Faro - Sacher         1.751 0.944 30   1.856  0.6966
 Faro - Slavon         1.751 1.090 30   1.607  0.8345
 Faro - Slovan        -0.857 0.944 30  -0.908  0.9950
 Faro - Tonkar         0.834 0.944 30   0.884  0.9959
 Flausist - Sacher     0.660 0.703 30   0.939  0.9936
 Flausist - Slavon     0.660 0.890 30   0.742  0.9989
 Flausist - Slovan    -1.948 0.703 30  -2.770  0.1924
 Flausist - Tonkar    -0.257 0.703 30  -0.365  1.0000
 Sacher - Slavon       0.000 0.944 30   0.000  1.0000
 Sacher - Slovan      -2.609 0.771 30  -3.385  0.0531
 Sacher - Tonkar      -0.917 0.771 30  -1.190  0.9683
 Slavon - Slovan      -2.609 0.944 30  -2.764  0.1944
 Slavon - Tonkar      -0.917 0.944 30  -0.972  0.9918
 Slovan - Tonkar       1.691 0.771 30   2.195  0.4814

Note: contrasts are still on the log(mu + 1) scale. Consider using
      regrid() if you want contrasts of back-transformed estimates. 
Degrees-of-freedom method: satterthwaite 
P value adjustment: tukey method for comparing a family of 10 estimates 

Effect sizes:

eff_size(emm,
         sigma = sigma(mod.spreizer),
         edf = df.residual(mod.spreizer))
Since 'object' is a list, we are using the contrasts already present.
 contrast             estimate    SE df lower.CL upper.CL
 (Astein - Bestein)      0.000 0.837 30  -1.7087   1.7087
 (Astein - Beton)       -0.148 0.517 30  -1.2036   0.9072
 (Astein - Capello)     -1.365 0.613 30  -2.6174  -0.1124
 (Astein - Faro)        -2.273 1.140 30  -4.5942   0.0490
 (Astein - Flausist)    -0.857 0.739 30  -2.3666   0.6528
 (Astein - Sacher)       0.000 0.854 30  -1.7450   1.7450
 (Astein - Slavon)       0.000 1.110 30  -2.2651   2.2651
 (Astein - Slovan)      -3.385 0.951 30  -5.3276  -1.4427
 (Astein - Tonkar)      -1.190 0.852 30  -2.9297   0.5489
 (Bestein - Beton)      -0.148 0.753 30  -1.6861   1.3897
 (Bestein - Capello)    -1.365 0.822 30  -3.0442   0.3144
 (Bestein - Faro)       -2.273 1.260 30  -4.8496   0.3044
 (Bestein - Flausist)   -0.857 0.920 30  -2.7359   1.0220
 (Bestein - Sacher)      0.000 1.010 30  -2.0728   2.0728
 (Bestein - Slavon)      0.000 1.240 30  -2.5262   2.5262
 (Bestein - Slovan)     -3.385 1.100 30  -5.6267  -1.1437
 (Bestein - Tonkar)     -1.190 1.010 30  -3.2584   0.8775
 (Beton - Capello)      -1.217 0.486 30  -2.2088  -0.2247
 (Beton - Faro)         -2.124 1.070 30  -4.3119   0.0631
 (Beton - Flausist)     -0.709 0.640 30  -2.0148   0.5973
 (Beton - Sacher)        0.148 0.753 30  -1.3897   1.6861
 (Beton - Slavon)        0.148 1.030 30  -1.9614   2.2578
 (Beton - Slovan)       -3.237 0.868 30  -5.0101  -1.4639
 (Beton - Tonkar)       -1.042 0.766 30  -2.6057   0.5212
 (Capello - Faro)       -0.908 1.080 30  -3.1050   1.2896
 (Capello - Flausist)    0.508 0.693 30  -0.9081   1.9241
 (Capello - Sacher)      1.365 0.822 30  -0.3144   3.0442
 (Capello - Slavon)      1.365 1.080 30  -0.8499   3.5798
 (Capello - Slovan)     -2.020 0.846 30  -3.7480  -0.2925
 (Capello - Tonkar)      0.174 0.802 30  -1.4637   1.8126
 (Faro - Flausist)       1.416 1.170 30  -0.9740   3.8054
 (Faro - Sacher)         2.273 1.260 30  -0.3044   4.8496
 (Faro - Slavon)         2.273 1.450 30  -0.6814   5.2267
 (Faro - Slovan)        -1.113 1.230 30  -3.6322   1.4071
 (Faro - Tonkar)         1.082 1.230 30  -1.4365   3.6008
 (Flausist - Sacher)     0.857 0.920 30  -1.0220   2.7359
 (Flausist - Slavon)     0.857 1.160 30  -1.5129   3.2267
 (Flausist - Slovan)    -2.528 0.973 30  -4.5162  -0.5403
 (Flausist - Tonkar)    -0.334 0.914 30  -2.2001   1.5330
 (Sacher - Slavon)       0.000 1.240 30  -2.5262   2.5262
 (Sacher - Slovan)      -3.385 1.100 30  -5.6267  -1.1437
 (Sacher - Tonkar)      -1.190 1.010 30  -3.2584   0.8775
 (Slavon - Slovan)      -3.385 1.310 30  -6.0516  -0.7188
 (Slavon - Tonkar)      -1.190 1.240 30  -3.7127   1.3318
 (Slovan - Tonkar)       2.195 1.040 30   0.0665   4.3230

sigma used for effect sizes: 0.7706 
Degrees-of-freedom method: inherited from satterthwaite when re-gridding 
Confidence level used: 0.95 

Plot

plot.spreizer <- dat.w  |> 
  # make plot
  mutate(jit = jitter(as.numeric(boar), 0.3)) |>  
  ggplot(aes(y   = prop.spreizer)) +
  geom_boxplot(aes(x   = boar, 
                   fill = husbandry),
               col = "black",
               outlier.shape = NA, width = 0.5) +
  geom_jitter(aes(x   = boar,
                  fill = husbandry),
              col = "black",
              position = position_jitterdodge(jitter.width  = 0.15,
                                              dodge.width   = 0.5),
              size = 2, alpha = 0.7,
              show.legend = FALSE) +
  scale_fill_manual(labels = c("Conventional",
                               "Ecological"),
                    values = c("#FFA040", "#008000")) +
  scale_y_continuous(lim = c(0, 20), 
                     breaks = seq(0, 20, 5),
                     minor_breaks = seq(0, 20, by = 1) ) +
  labs(x = "",
       y = "Proportion spreizer [%]",
       fill = "Husbandry") +
  my_theme +
  theme(legend.position = "top") +
  guides(y = guide_axis(minor.ticks = TRUE),
         fill = guide_legend(override.aes = list(linetype = 0, size  = 5)))
plot.spreizer

Proportion low body weight

Model

hist(dat.w$prop.lbw,
     breaks = 30)

hist(log(dat.w$prop.lbw+1),
     breaks = 30)

mod.lbw <- lmerTest::lmer(log(dat.w$prop.lbw+1) ~ 
                                   boar +
                                   # random intercept for sows
                                   (1 | sow),
                                 data    = dat.w,
                                 REML    = TRUE,
                                 control = contr)
summary(mod.lbw)
Linear mixed model fit by REML. t-tests use Satterthwaite's method [
lmerModLmerTest]
Formula: log(dat.w$prop.lbw + 1) ~ boar + (1 | sow)
   Data: dat.w
Control: contr

REML criterion at convergence: 112.5

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-1.86969 -0.18804  0.02804  0.51646  1.57564 

Random effects:
 Groups   Name        Variance Std.Dev.
 sow      (Intercept) 0.03886  0.1971  
 Residual             1.73722  1.3180  
Number of obs: 40, groups:  sow, 15

Fixed effects:
             Estimate Std. Error      df t value Pr(>|t|)    
(Intercept)    3.2141     0.5959 29.9846   5.394 7.67e-06 ***
boarBestein   -1.0909     1.1147 29.9854  -0.979    0.336    
boarBeton     -0.7084     0.6860 28.0416  -1.033    0.311    
boarCapello    0.2217     0.7779 28.1133   0.285    0.778    
boarFaro       0.6673     1.4593 29.9971   0.457    0.651    
boarFlausist  -1.1866     0.9703 28.3199  -1.223    0.231    
boarSacher    -1.2459     1.1145 29.9964  -1.118    0.272    
boarSlavon    -0.3499     1.4594 29.9994  -0.240    0.812    
boarSlovan     0.4711     1.1147 29.9849   0.423    0.676    
boarTonkar    -0.2093     1.1078 22.6429  -0.189    0.852    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) brBstn borBtn brCpll boarFr brFlss brSchr boarSlavn
boarBestein -0.534                                                    
boarBeton   -0.863  0.465                                             
boarCapello -0.759  0.406  0.658                                      
boarFaro    -0.408  0.218  0.358  0.310                               
boarFlausst -0.610  0.326  0.529  0.465  0.249                        
boarSacher  -0.534  0.286  0.465  0.412  0.218  0.326                 
boarSlavon  -0.408  0.218  0.355  0.310  0.167  0.249  0.218          
boarSlovan  -0.534  0.286  0.463  0.406  0.218  0.332  0.286  0.230   
boarTonkar  -0.526  0.281  0.458  0.399  0.215  0.321  0.281  0.215   
            boarSlovn
boarBestein          
boarBeton            
boarCapello          
boarFaro             
boarFlausst          
boarSacher           
boarSlavon           
boarSlovan           
boarTonkar   0.281   
round(drop1(mod.lbw, test = 'Chisq'), 3)
Single term deletions using Satterthwaite's method:

Model:
log(dat.w$prop.lbw + 1) ~ boar + (1 | sow)
     Sum Sq Mean Sq NumDF  DenDF F value Pr(>F)
boar 11.706   1.301     9 28.663   0.749  0.662
car::Anova(mod.lbw,
           test.statistic = "Chisq",
           type = 2)
Analysis of Deviance Table (Type II Wald chisquare tests)

Response: log(dat.w$prop.lbw + 1)
      Chisq Df Pr(>Chisq)
boar 6.7386  9     0.6643
car::Anova(mod.lbw,
           test.statistic = "F",
           type = 2)
Analysis of Deviance Table (Type II Wald F tests with Kenward-Roger df)

Response: log(dat.w$prop.lbw + 1)
          F Df Df.res Pr(>F)
boar 0.6549  9 27.009 0.7409

Model diagnostics

performance::check_model(mod.lbw)

Plot

plot.lbw <- dat.w  |> 
  # make plot
  mutate(jit = jitter(as.numeric(boar), 0.3)) |>  
  ggplot(aes(y   = prop.lbw)) +
  geom_boxplot(aes(x   = boar, 
                   fill = husbandry),
               col = "black",
               outlier.shape = NA, width = 0.5) +
  geom_jitter(aes(x   = boar,
                  fill = husbandry),
              col = "black",
              position = position_jitterdodge(jitter.width  = 0.15,
                                              dodge.width   = 0.5),
              size = 2, alpha = 0.7,
              show.legend = FALSE) +
  scale_fill_manual(labels = c("Conventional",
                               "Ecological"),
                    values = c("#FFA040", "#008000")) +
  scale_y_continuous(lim = c(0, 82), 
                     breaks = seq(0, 82, 20),
                     minor_breaks = seq(0, 82, by = 5) ) +
  labs(x = "",
       y = "Proportion LBW [%]",
       fill = "Husbandry") +
  my_theme +
  theme(legend.position = "top") +
  guides(y = guide_axis(minor.ticks = TRUE),
         fill = guide_legend(override.aes = list(linetype = 0, size  = 5)))
plot.lbw

Combined plots: Figure

# Combine plots with a designated area for the legend
combined <- (plot.bodyweight + 
             plot.piglets + 
             guide_area() + 
             plot.males + 
             plot.stillborn+
             plot.spreizer) + 
  plot_layout(guides = "collect") +
  plot_annotation(tag_levels = 'A') &
  theme(plot.tag = element_text(face = "bold", size = 20),
        legend.position = "right",
        legend.direction = "vertical")
combined

png("./plots/figure-boar.png",
     width = 400, height = 300, units = "mm",
     pointsize = 10, res = 600)

combined

dev.off()
png 
  2 

How to cite R

“All analyses were performed using R Statistical Software (version 4.4.2; R Core Team 2024)”.

Reference: R Core Team (2024). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. URL https://www.R-project.org/.

citation()
To cite R in publications use:

  R Core Team (2024). _R: A Language and Environment for Statistical
  Computing_. R Foundation for Statistical Computing, Vienna, Austria.
  <https://www.R-project.org/>.

Ein BibTeX-Eintrag für LaTeX-Benutzer ist

  @Manual{,
    title = {R: A Language and Environment for Statistical Computing},
    author = {{R Core Team}},
    organization = {R Foundation for Statistical Computing},
    address = {Vienna, Austria},
    year = {2024},
    url = {https://www.R-project.org/},
  }

We have invested a lot of time and effort in creating R, please cite it
when using it for data analysis. See also 'citation("pkgname")' for
citing R packages.
version$version.string
[1] "R version 4.4.2 (2024-10-31 ucrt)"
citation("tidyverse")
Um Paket 'tidyverse' in Publikationen zu zitieren, nutzen Sie bitte:

  Wickham H, Averick M, Bryan J, Chang W, McGowan LD, François R,
  Grolemund G, Hayes A, Henry L, Hester J, Kuhn M, Pedersen TL, Miller
  E, Bache SM, Müller K, Ooms J, Robinson D, Seidel DP, Spinu V,
  Takahashi K, Vaughan D, Wilke C, Woo K, Yutani H (2019). "Welcome to
  the tidyverse." _Journal of Open Source Software_, *4*(43), 1686.
  doi:10.21105/joss.01686 <https://doi.org/10.21105/joss.01686>.

Ein BibTeX-Eintrag für LaTeX-Benutzer ist

  @Article{,
    title = {Welcome to the {tidyverse}},
    author = {Hadley Wickham and Mara Averick and Jennifer Bryan and Winston Chang and Lucy D'Agostino McGowan and Romain François and Garrett Grolemund and Alex Hayes and Lionel Henry and Jim Hester and Max Kuhn and Thomas Lin Pedersen and Evan Miller and Stephan Milton Bache and Kirill Müller and Jeroen Ooms and David Robinson and Dana Paige Seidel and Vitalie Spinu and Kohske Takahashi and Davis Vaughan and Claus Wilke and Kara Woo and Hiroaki Yutani},
    year = {2019},
    journal = {Journal of Open Source Software},
    volume = {4},
    number = {43},
    pages = {1686},
    doi = {10.21105/joss.01686},
  }
citation("patchwork")
Um Paket 'patchwork' in Publikationen zu zitieren, nutzen Sie bitte:

  Pedersen T (2024). _patchwork: The Composer of Plots_. R package
  version 1.3.0, <https://CRAN.R-project.org/package=patchwork>.

Ein BibTeX-Eintrag für LaTeX-Benutzer ist

  @Manual{,
    title = {patchwork: The Composer of Plots},
    author = {Thomas Lin Pedersen},
    year = {2024},
    note = {R package version 1.3.0},
    url = {https://CRAN.R-project.org/package=patchwork},
  }
citation("lmerTest")
To cite lmerTest in publications use:

  Kuznetsova A, Brockhoff PB, Christensen RHB (2017). "lmerTest
  Package: Tests in Linear Mixed Effects Models." _Journal of
  Statistical Software_, *82*(13), 1-26. doi:10.18637/jss.v082.i13
  <https://doi.org/10.18637/jss.v082.i13>.

Ein BibTeX-Eintrag für LaTeX-Benutzer ist

  @Article{,
    title = {{lmerTest} Package: Tests in Linear Mixed Effects Models},
    author = {Alexandra Kuznetsova and Per B. Brockhoff and Rune H. B. Christensen},
    journal = {Journal of Statistical Software},
    year = {2017},
    volume = {82},
    number = {13},
    pages = {1--26},
    doi = {10.18637/jss.v082.i13},
  }
citation("car")
To cite the car package in publications use:

  Fox J, Weisberg S (2019). _An R Companion to Applied Regression_,
  Third edition. Sage, Thousand Oaks CA.
  <https://www.john-fox.ca/Companion/>.

Ein BibTeX-Eintrag für LaTeX-Benutzer ist

  @Book{,
    title = {An {R} Companion to Applied Regression},
    edition = {Third},
    author = {John Fox and Sanford Weisberg},
    year = {2019},
    publisher = {Sage},
    address = {Thousand Oaks {CA}},
    url = {https://www.john-fox.ca/Companion/},
  }
citation("emmeans")
Um Paket 'emmeans' in Publikationen zu zitieren, nutzen Sie bitte:

  Lenth R (2024). _emmeans: Estimated Marginal Means, aka Least-Squares
  Means_. R package version 1.10.6,
  <https://CRAN.R-project.org/package=emmeans>.

Ein BibTeX-Eintrag für LaTeX-Benutzer ist

  @Manual{,
    title = {emmeans: Estimated Marginal Means, aka Least-Squares Means},
    author = {Russell V. Lenth},
    year = {2024},
    note = {R package version 1.10.6},
    url = {https://CRAN.R-project.org/package=emmeans},
  }
citation("performance")
Um Paket 'performance' in Publikationen zu zitieren, nutzen Sie bitte:

  Lüdecke et al., (2021). performance: An R Package for Assessment,
  Comparison and Testing of Statistical Models. Journal of Open Source
  Software, 6(60), 3139. https://doi.org/10.21105/joss.03139

Ein BibTeX-Eintrag für LaTeX-Benutzer ist

  @Article{,
    title = {{performance}: An {R} Package for Assessment, Comparison and Testing of Statistical Models},
    author = {Daniel Lüdecke and Mattan S. Ben-Shachar and Indrajeet Patil and Philip Waggoner and Dominique Makowski},
    year = {2021},
    journal = {Journal of Open Source Software},
    volume = {6},
    number = {60},
    pages = {3139},
    doi = {10.21105/joss.03139},
  }

Session Info

sessionInfo()
R version 4.4.2 (2024-10-31 ucrt)
Platform: x86_64-w64-mingw32/x64
Running under: Windows 11 x64 (build 26100)

Matrix products: default


locale:
[1] LC_COLLATE=German_Germany.utf8  LC_CTYPE=German_Germany.utf8   
[3] LC_MONETARY=German_Germany.utf8 LC_NUMERIC=C                   
[5] LC_TIME=German_Germany.utf8    

time zone: Europe/Berlin
tzcode source: internal

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] performance_0.13.0 emmeans_1.10.6     car_3.1-3          carData_3.0-5     
 [5] lmerTest_3.1-3     lme4_1.1-36        Matrix_1.7-1       patchwork_1.3.0   
 [9] lubridate_1.9.4    forcats_1.0.0      stringr_1.5.1      dplyr_1.1.4       
[13] purrr_1.0.2        readr_2.1.5        tidyr_1.3.1        tibble_3.2.1      
[17] ggplot2_3.5.1      tidyverse_2.0.0   

loaded via a namespace (and not attached):
 [1] gtable_0.3.6        bayestestR_0.15.1   xfun_0.50          
 [4] htmlwidgets_1.6.4   ggrepel_0.9.6       insight_1.0.1      
 [7] lattice_0.22-6      tzdb_0.4.0          numDeriv_2016.8-1.1
[10] vctrs_0.6.5         tools_4.4.2         Rdpack_2.6.2       
[13] generics_0.1.3      datawizard_1.0.0    parallel_4.4.2     
[16] pbkrtest_0.5.3      sandwich_3.1-1      pkgconfig_2.0.3    
[19] lifecycle_1.0.4     compiler_4.4.2      farver_2.1.2       
[22] munsell_0.5.1       codetools_0.2-20    htmltools_0.5.8.1  
[25] yaml_2.3.10         Formula_1.2-5       pillar_1.10.1      
[28] nloptr_2.1.1        MASS_7.3-61         reformulas_0.4.0   
[31] boot_1.3-31         abind_1.4-8         multcomp_1.4-26    
[34] nlme_3.1-166        tidyselect_1.2.1    digest_0.6.37      
[37] mvtnorm_1.3-3       stringi_1.8.4       labeling_0.4.3     
[40] splines_4.4.2       fastmap_1.2.0       grid_4.4.2         
[43] colorspace_2.1-1    cli_3.6.3           magrittr_2.0.3     
[46] survival_3.7-0      broom_1.0.7         TH.data_1.1-3      
[49] withr_3.0.2         backports_1.5.0     scales_1.3.0       
[52] timechange_0.3.0    estimability_1.5.1  rmarkdown_2.29     
[55] zoo_1.8-12          hms_1.1.3           coda_0.19-4.1      
[58] evaluate_1.0.3      knitr_1.49          rbibutils_2.3      
[61] mgcv_1.9-1          rlang_1.1.4         Rcpp_1.0.13-1      
[64] xtable_1.8-4        glue_1.8.0          see_0.9.0          
[67] rstudioapi_0.17.1   minqa_1.2.8         jsonlite_1.8.9     
[70] R6_2.6.1